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Abstract 

In the current work, first of all, a general theoretical two- 
phase model was proposed to simulate the flow in the resin 
transfer molding RTM. This model has been widely used in 
geologic media and petroleums reservoirs and satisfied 
results have been given. The similarity between these 
medium and our composite domain either in heterogeneity 
or in fluid flow has motivated this choice. Based on 
multiphase darcy's law, the model has been developed to 
simulate the saturation distribution for two compressible 
immiscible phases which are resin and the air, it may control 
the formation of the void in composite material. Then, a 
particular case of the model was simulated numerically 
which is well known as the buckely-levrett model. Using the 
first time the control volume -finite element CV/FEM method, 
the results of this study agree qualitatively with 
e xpe rime ntal and nume rical findings . 
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Introduction 

Resin Transfer Molding (RTM) is a closed molding 
process in which a thermodurcissable resin is injected 
into a preplaced fiber preform. RTM offers the 
promise of producing low cost composite parts with 
complex structures and large near net shapes. 
Relatively fast cycle times with good surface definition 
and appearance are easily achievable. For these 
reasons, this process has gained the attention of 
composite manufacturers especially in automotive 
industries (Owen et al., 1991), and the subject of 
numerous experimental and numerical studies has 
been done (Hattabi et al., 2005, 2008; Saouab et al.. 


2006; Samir et al. 2009, 2011). Nevertheless, almost all 
of these numerical studies have considered the 
problem as a one phase flow model which takes into 
account only the resin flow, therefore all parameters 
determined rely on this phase. While in reality, RTM is 
a technique where the motion of two phases plays 
significant role and affects the process parameters. The 
motivation for the research work presented in this 
paper arises from the need to develop an efficient and 
accurate model for RTM. This model accounts for the 
multiphase nature and more other singularities of this 
process, and in the same time the presence and the 
interaction of the two different phases are taken into 
consideration; liquid (resin) and gaz (air), which can 
also be useful to treat the formation and migration of 
the void that affects the quality and the performance 
of the composite material. The two-phase general 
model developed in our work is extrapolated from soil 
mechanics (Bastian, 1999; Chen et al., 1994, 2006; 
Douglas et al., 1987; Geiger et al., 2003; Webb, 2006), 
which enable us to treat the case where one or the two 
phases are compressible, the permeability dependence 
on saturation, and even, the fiber porosity time- 
dependent. In spite of the obvious advantages offered 
by the two- phase model, few researches (Chui et al., 
1997; Pillai et al., 1996) have been performed in this 
direction. 

Aiming to test the validity of our model, a particular 
case of the filling process was simulated numerically, 
based on the assumption that both of the two phases 
are incompressible, which is well known as the 
Buckely-levrett model (Chui et al. 1997), the originality 
of this study was from the implementation of the 
CV/FEM methods to track the flow in the molding. A 
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node-centred control volume coupled with a finite 
element method has been depicted on an unstructured 
triangular grid which can be used to accurately and 
efficiently model multi-phase flow in porous media 
(Timoorei, 2005; Chen, 2006). The results of our 
simulation are qualitatively quite the same with 
experimental (Hattabi et al., 2005, 2008; Saouab et al., 
2006) and numerical findings (Samir et al., 2009,2011). 


Mathematical Model 


The equations describing the flow of two compressible 
and immiscible fluids can be derived by combining 
Darcy's law and mass conservation equation written 
for each phase individually. The phase flow equations 
are given in a fractional flow formulation, i.e, in terms 
of saturation and a global pressure. The main reason 
for this fractional flow approach is that efficient 
numerical methods can be devised to take advantage 
of many physical properties inherent to the flow 
equations. 


The parabolic equation describing the fluid pressure/? 
in the mold is then defined as: 


(c a S a + c r S r ) J ^ + q be, S r ) 

dS 
dt 
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In addition to the above system, a law of state is added 
to describe the behaviour of each phase: 

S r + S a = 1 (3) 

Where 0 < S r ,S a < 1. 


Equation (4) signifies that the two phases fill the entire 
pore of the porous media. 


In these equations, subscript a and r refer respectively 
to air and resin, S is the saturation, K is the absolute 
permeability of the fiber mat, 0 is the porosity of the 
media, and c is the compressibility of each phase, A is 
the mobility, p is the density and g is the acceleration 
due to gravity. 


With:- = 4> — + — V. 

dt dt Sn 
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In equations (1) and (2), the average fluid pressure p is 
given by: 

P = "(Pr + Pa> (5) 

Vt is the total velocity: 


K = K + K ( 6 ) 

The difference of pressures between the two phases. 


named the capillary pressure p c is given by: 

Pc = Pa Pr ( 7 ) 

A t is the total mobility: 

i t t k rr (s r ) ( k r a ( S a ) /ox 

A t - A r + A a - — + " (o) 

Mr Ma 

/c r r (5 r ), k r a (5 a ) denote relative permeabilities which 
take values between 0 and 1, strongly depending on 
the saturation Sr and Sa. This dependence is not yet 
well established for general porous media, and the 
experimental research is therefore necessary to 
develop this relationship. 

fr represents the flow fraction function of the resin: 

fc r , r (5 r ) 

fr = a 7 = fc r.r(SrAr.a(^T (9) 

Mr Ma 

In our two-phase model, it is considered that no phase 
transfer between the two fluids will occur. It's also 
assumed that the capillary and the gravitational force 
can be neglected. Thus, the flow equations become: 


o 

II 

> 

(10) 

K = -K A t Vp 

(11) 

4>^+v(f r v t ) =0 

(12) 


This system is well known by the Buckley -levrett 
equation (Chui et al., 1997). It is more used in 
modelling the geologic media (Chen et al., 1994, Chui 
et al., 1997), and recently in few works in composite 
field (Chui et al., 1997). 

The initial conditions are resin and air saturations. For 
the boundary conditions, a Dirichlet boundary 
condition of pressure at the inlet and the outlet gate: 

pinlet = pinj (13) 

pfront = 0 (- 1 ^) 

Whereas a Neumann boundary condition (no flow) is 
assumed at the other boundary: 

v.n - 0 (15) 

Where v is the velocity and n is the normal to the mold 
wall. 

The relative permeability dependence of the saturation 
can be modelled as (Chui et al., 1997): 

K,a^J=S 2 a . (16) 

k rr (S r ) = S 2 . (17) 


Numerical Procedure 

Since the system of equation of the two-phase model is 
rearranged in such a way to have a form similar to that 
of the one phase model equation (Samir et al. 2009, 
2011). To solve the Buckley -levrett system, many 
researches have been made based on level set method 
(Soukane et al., 2006), volume of fluid VOF, finite 
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difference and standard finite element method 
(Bruschke et al., Lin et.). In our study and contrarily to 
all methods used previously, the resin flow problem 
has been resolved by means of the control 
volume/finite element method (CV/FEM). This method 
was selected because of the advantage that this fixed 
mesh method eliminates the need for remeshing the 
resin-filled domain for each time step, thus the flow 
simulation for a complex geometry can be done 
rapidly and efficiently. Numerous studies have shown 
that CV/FEM simulations (chen et al., 1994; Shojaeia et 
al., 2003, 2004, 2006) yield very good results in 
comparison to traditional finite difference (Trochu et 
al., 1992; Hattabi et al., 2005, 2008) or finite element 
methods (Chang et al., 1998, 2003). 

The use of the CV/FEM has many advantages. In 
particular: 

- The geometric flexibility of the finite-element method 
allows for large variations in scale to be modelled 
efficiently. 

- The parabolic equation for the fluid pressure and the 
hyperbolic conservation equation can be solved 
efficiently by the finite element method. 

In addition to the geometric flexibility of using 
unstructured meshes, this method has the advantage 
that only the fluid pressure field needs to be computed 
by the finite element method. 

Finit Element Method 

The finite element method is used to approximate the 
spatial operators. This yields the fluid pressure at the 
nodes of the finite elements at each time step. The 
pressure is calculated using the Galerkin approach, 
and the governing equation of the flow model is: 

V(/a t Vp) = 0 (18) 

Using the procedure of Galerkin to this equation, we 
obtain: 

4 V (O t Vp)dfl = 0 (19) 

J n KAtVpVQidn = 0 (20) 

For all basis functions </>;, the pressure p has the 
decomposition: 

p(x, t) = Y?- =1 Pj (t) 0 ; (x) (21) 

The equations can then be written as follows: 

M{Pj)= 0 ( 22 ) 

With: 

K ij = 4 XMMjdCl = 0 (23) 

After the distribution of the pressure is determined in 


the domain, the nodal fluid pressures can be 
differentiated to compute the total velocities at the 
centre of the finite elements, using Darcy's law 
equation (eqt. 11). 

Control Volume Technique 

The control volume technique is a numerical approach 
to track flow front location in fluid dynamics problems 
involving flow with a free-surface (Chen et al., 2006) 

In the control volume approach, the molding is first 
discretized into finite elements. By subdividing the 
elements into smaller sub- volumes, a control volume 
is constructed around each node (figure 1). The 
concept of fill factor is introduced to monitor the fluid 
volume in each control volume; and defined as the 
ratio of the fluid volume in the control volume to the 
total volume of the control volume. The fill factor 
takes values from 0 to 1, of which 0 represents that the 
control volume is empty and 1 when it is filled. The 
control volumes can be empty, partially filled, or 
completely filled. The numerical flow front is 
construction of the nodes that have partially filled 
control volumes. At each time step, fill factors are 
calculated based on the resin velocity and flow into 
each nodal control volume. If the resin does not gel, 
the flow front is thus updated until all the control 
volumes are full and the infiltration is complete. 



FIG. 1 ONE CONTROL VOLUME Vi AT THE NODE z, WITH 
TRIANGULAR ELEMENT ai -eis 

There are two basic components to the control volume 
technique: (1) control volume formation, (2) flow front 
advancement. 

1) Control Volume Formation 

The molding cavity is first discretized using finite 
elements, and then each individual element is 
further divided into sub- volumes; each of whcih is 
associated with one of the nodes on the element 
(figure 1). The number of the sub -volumes is equal 
to that of nodes in the element. The control volume 
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for a particular node is composed of all of the sub- 
volumes associated with that node. 


2) Flow Front Advancement 


The concept of fill factor introduced in one-phase 
model to monitor the progression of the flow front 
is replaced in our two-phase model by the 
saturation Sr of the wetting phase (resin). 


The equation that leads to updating the saturation 
Sr is the same as it describes the mass conservation 
of the resin: 

tpjf+VW = 0 (24) 

After calculating the total velocity, this one is 
introduced in the equation (24) to determine the 
saturation for the next iteration. The integration of 
equation (24) over an arbitrary volume Vi and the 
application of the divergence theorem yields: 

4>$v^ dV i = - f v Xf r V t )dV t (25) 

Within each control volume Vi, Sr is constant. 
Discretization of Equation (25) using Euler's 
method leads to: 

5 "‘ + 1 = 3 * S "i + P V tj n J ) (26) 

k 

Where Z is the summation over all segments j 

belonging to the control volume. At is the time step, 
Ai is the area of the control volume and nj is the 
outward normal vector to segment j. 


The saturation Sr takes values from 0 to 1, where 0 
represents totally empty and 1 totally full. When 
the control volume is partially full, the saturation 
Sr varies between 0 and 1. The numerical flow front 
is construction of the nodes that have partially full 
control volumes as shown in fig. 2. 



Empty CV 





Numerical 
flow front 


FIG. 2 NUMERICAL AND ACTUAL FLOW FRONT 


3) Time Step Calculation 

The time step for the next iteration must be 

calculated before the solution can proceed. The 


optimal time step would be when the fluid just fills 
one control volume. If a larger step was chosen, the 
flow front would over -run the control volume and 
result in a loss of mass from the system. 


The flow charts of the numerical schemes 
developed in this study are illustrated in fig. 3. 



FIG. 3 FLOW CHART OF THE NUMERICAL PROCEDURE OF 
FE/CV TECHNIQUE 


Results and discussion 

In the simulation of the RTM filling process, the resin 
is injected into a square molding cavity (4000 x 4000) 
mm2, and the inlet condition can be either a constant 
volumetric flow rate or a constant injection pressure 
for the resin. Here constant pressure is taken as the 
inlet boundary condition. The details of the numerical 
procedure are described by the flow chart in figure 3. 
The physical properties of the fiber mats and resin are 
listed in table 1. 

Table 1 Physical properties of fiber and processing conditions 

USED IN NUMERICAL STUDY 


RTM parameters 

Permeability 

Kxx = IO 9 Kyy=l(P 

Resin viscosity 

jur =0.05 Pa.s 

Air viscosity 

|Lla =1.8.10' 5 Pa.s 

Inj e c tio n pre s sure 

pnj. = 4 x 10 5 Pa 

Porosity 

O =0,696 


When the molding is completely filled, the saturation 
Sr is the same in all nodes and equal to 1. Then the 
equation (18) is reduced to the one considered in one- 
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phase flow (Hattabi et al., 2005, 2008; Samir et al., 2009, 
2011), therefore it has the same solution given in (Cai) 
by: 

Where P0 is the pressure at the inlet, x is the flow front 
position. 

Figure 4 shows a comparison between our numerical 
result and the analytical result of pressure presented 
in equation (28). An excellent agreement is observed. 



FIG. 4 PRESSURE VARIATION OVER THE MOLD LENGTH 

Figure 5 is a classical observation that comes to 
reinforcing the validity of the present model, in 
addition, it has been demonstrated that the viscosity is 
among other parameters that affect the cycle time of 
the process. Indeed, the time of filling increases with 
high resin viscosity. So, in order to optimize the 
production cost, it's recommended to minimize the 
time of the cycle by using a not very viscous resin. 



FIG. 5 EVOLUTION OF MOLD TIME FILLING TIME WITH RESIN 
VISCOSITY 

The front position presented in figure 6 has a 
quadratic pattern in the course of the time, witch 
confirms for another time the conformity of the model 
with the previous findings (Saad et al., 2011, 2012; 


Samir et al., 2009,2011). 



In this model, the total time of the filing becomes 
higher 895.39 s as compared with the model of one 
phase flow 500.58 s (figure 7), because in our model 
both the motion of resin and air are taken into account, 
and this increment is principally due to the resistance 
proved by air to the resin advancement. 



FIG. 7 FRONT POSITION COMPARISON BETWEEN TWO 
MODELS 



FIG. 8 SCHEMATIC OF THE COMPLEX MOULD GEOMETRY 
WITH SINGLE GATE AND BOUNDARY CONDITIONS 

The modelling of complex geometry is a more realistic 
concern in RTM application field, that's why it is 
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intended to simulate with our model the filling of a 
complex shape, which has the specificity including an 
insert inside the piece. The filling of the mould has 
been done under constant injection pressure with only 
a single gate injection, and five events where the 
pressure is equal to zero (fig. 8). 

The distribution of pressure in the course of molding 
filling is also affected in the current model (figure 9). 
This confirms again the importance of taking into 
consideration the motion of the two phases (resin and 
air) in order to define more precisely all process 
parameters that allow a final product with high 
performance and low cost. 

Pressure 

FI 4 

— 3.75 



FIG. 9 NUMERICAL PRESSURE DISTRIBUTION IN (a) ONE 
PHASE FLOW (b) TWO PHASES FLOW 



FIG. 10 RESIN SATURATION PROFILE AT DIFFERENT TIME 

Figure 10 shows the saturation profile as a function of 
the distance from inlet, at various time of filling, and 
this result agree well with the saturation pattern of 


Chang et al. At the initial time, the saturation is equal 
to one just at the inlet gate and zero at the flow front. 
The flow field expands gradually during the filling 
process, and the wetting phase saturation Sr increases 
gradually from 0 to the final steady state value 1, to 
completely fill the entire molding. This implies that 
the voids will disappear if the resin is continuously 
filled. 

Conclusion and outlook 

The theoretical general model seems to be promising 
in the simulation of resin transfer molding, since it 
takes into consideration almost all particularity of 
RTM filling process. This constatation comes from the 
fact that the Buckley- lev rett model, a particular case of 
this model, gave a satisfied result in the case of two 
incompressible steady states problem. 

The main result of this work is that the Buckley -lev rett 
model presented in a form similar to that of a one 
phase model, can be resolved by the CV/FEM method 
to predict accurately all process variables. 
Furthermore, the control volume finite element 
method applied for the first time to a two-phase flow 
model allows one to predict with a cost-effective tool 
and a high accuracy the pressure evolution over the 
saturated region of the mould at each time step and to 
track the flow front position during molding filling, 
and the pressure can be also investigated in molding 
with simple and complex geometries including inserts. 

To improve the usefulness of this model, it is 
suggested in the next work, to simulate numerically 
this theoretical model, by taking into account the 
compressibility of the two phases. 

REFERENCES 

Bastian, Peter. Numerical Computation of Multiphase Flows 
in Porous Media, doctorat d'habilite (1999):16-27. 
Bruschke, M.V. and Advani, S.G. A finite element/control 
volume approach to mold filling in anisotropic porous 
media. Polymer Composites., 11 (1990): 398-405. 

Cai, Z. Simplified mold filling simulation in resin transfer 
molding. Journal of Composite Materials, 26 (1992): 
2606-2629. 

Chang, C.Y and Shih, M. S. Numerical simulation on the 
void distribution in the fiber mats during the filling stage 
of RTM, Journal of reinforced plastics and composite, 22 
(16) (2003): 1437-1454. 


63 



www. s e ipub . o r g/r pim 


Re search on Precision Instrument and Machinery Vol. 2 Iss. 4, December 2013 


Chang, C.Y. and Hourng, L.W. Numerical simulation for the 
transverse impregnation in resin transfer moulding. 
Journal of reinforced plastics and composite, 17 (2) (1998): 
165-182. 

Chen, Z., Es pedal, M. and Ewing. R. E. Finite element 
analysis of multiphase flow in groundwater hydrology. 
Institute of mechanic and its applications preprint series 
# 1268 University of Minnesota (1994). 

Chen, Zhangxin. On the control volume finite element 
methods and their applications to multiphase flow, 
networks and heterogeneous media, 1 (4) (2006): 689-706. 

Chui, W. Ketal., Process modeling in resin transfer molding 
as a method to enhance product quality, SIAM review., 
39 (4) (1997): 714-727. 

Douglas, J. and Yirang, Y. Numerical simulation of 
immiscible flow in porous media based on combining the 
method of characteristic with mixed finite element 
procedure. Institute of mechanic and its applications 
pre print s e rie s # 318 Unive rsity of Minne s ota (1987). 

Geiger, S., Robertsy, S., Matthaiz, S. K. and Zoppoux, C. 
Combining Finite Volume and Finite Element Methods 
to Simulate Fluid Flow in Geologic Media, ANZIAM J., 
44 (E) (2003): C180-C201. 

Hattabi, M., Echaabi, J.and Bensalah, M. O. Numerical and 
experimental analysis of the resin transfer molding 
process, Korea -Australia Rheology Journal., 20 (1) (2008): 
7-14. 

Hattabi, M., Echaabi, J., Bensaleh, O., Saouab, A. and Breard, 
J. Flow Analyzes during on-line and Radial Injections 
Application in the Measure of Permeability, Journal of 
Reinforced Plastics and Composites., 24 (2005): 1909-1920. 

Hattabi, M., Snaike, I., Echaabi, J. and Bensalah, M. O. 
Simulation du front d'ecoulement dans les procedes de 
moulage des composites liquides, Comptes rendus 
mecanique., 333 (7) (2005): 585-591. 

Lin, M., Hahn, H.T. and Huh, H. A finite element simulation 
of resin transfer molding based on partial nodal 
saturation and implicit time integration. Composites Part 
A, 29 (1998): 541-550. 

Owen, M.J, Middleton, V. and Hutcheon, K.F. The 
development of resin transfer molding (RTM) for volume 
manufacture, Proc I Mech E Design in Composite 
Materials (1989): 107-114. 


Owen, M.J, Rudd, C.D, Middleton, V., Kendall, K.N. and 
Revill, I.D. Resin transfer molding (RTM) for automotive 
components. Composite Material Technology ASME., 37 
(1991): 177-183. 

Pillai, KM. and Advani, S.G. Modeling of Void Migration in 
resin Transfer Molding Process, In: Proceedings of the 
1996 ASME International Mechanical Engineering 
Congress and Exhibition (Nov 17-22), Atlanta GA (1996) 
page 4. 

Saad, A., Echchelh, A., Hattabi, M. And El Ganaoui, M. A 
fast computational model to the simulation of non- 
isothermal mold filling process in resin transfer molding. 
Journal of polymer composites., 32 (6) (2011): 857-868. 

Saad, A., Echchelh, A., Hattabi, M. and El Ganaoui, M. An 
Improved Computational Method For Non Isothermal 
Resin Transfer Molding Simulation, Thermalsciences, 15, 
Issue suppl. 2 (2011): 275-289. 

Saad, A., Echchelh, A., Hattabi, M. and El Ganaoui, M. 
Numerical Simulation of Thickness Variation Effect on 
Resin Transfer Molding Process, Journal of Polymer 
composites., 33(1) (2012): 10-21. 

Samir, J., Hattabi, M. and Echaabi, J. Numerical Algorithm 
and Adaptive Meshing for Simulation the Effect of 
Variation Thickness in Resin Transfer Molding Process, 
Composites Part B., 42 (5) (2011): 1015-1028. 

Samir, J., Hattabi, M., Echaabi, J., Saouab, A. and Park, C. H 
Simulation of mold filling in RTM process by the control 
volume /finite element method, Africa in revue in the 
research of informatics and applied mathematics., 10 
(2009) : 1-15. 

Saouab, A. and Breard, J. Analytical modeling of CRTM and 
RTM Processes Part A: Mathematical developments. 
International Journal of Forming Processes., 9(2006): 309- 
344. 

Shojaeia A., Numerical simulation of three-dimensional flow 
and analysis of filling process in compression resin 
transfer moulding. Composites: Part A., 37 (2006): 1434- 
1450. 

Shojaeia, A., A numerical study of filling process through 
multilayer preforms in resin injection/compression 
molding. Composites Science and Technology., 66 (2006): 
1546-1557. 

Shojaeia, A., Ghaffariana, S.R. and Karimianb, S.M.H 


64 



Re s e arc h o n Pr e cis io n I ns tr ument and Ma chine ry V o 1. 2 1 s s . 4, Dec embe r 20 1 3 


www. s e ipub .or g/ r pim 


Numerical Analysis of Controlled Injection Strategies in 
Resin Transfer Molding, Journal of reinforced plastics 
and composite ., 37 (11) (2003): 1011-1035. 

Shojaeia, A., Ghaffariana, S.R. and Karimianb, S.M.H 

Numerical simulation of three-dimensional mold filling 
process in resin transfer molding using quasi-steady 
state and partial saturation formulations. Composites 
Science and Technology., 62 (2002): 861-879. 

Shojaeia, A., Ghaffariana, S.R. and Karimianb, S.M.H. 

Numerical Simulation of Three-Dimensional Mold 
Filling in Resin Transfer Molding, Journal of reinforced 
plastics and composite., 22 (16) (2003): 1497-1529. 

Shojaeia, A., Ghaffariana, S.R. and Karimianb, S.M.H. 

Simulation of the three-dimensional non-isothermal 
mold filling process in resin transfer molding. 

Composites Science and Technology., 63 (2003): 1931- 
1948. 

Shojaeia, A., Ghaffariana, S.R. and Karimianb, S.M.H Three- 


dimensional process cycle simulation of composite parts 
manufactured by resin transfer molding. Composite 
Structures., 65 (2004): 381-390. 

Soukane, S. and Trochu, F. Application of the level set 
method to the simulation of resin transfer molding. 
Composites Science and Technology, 66 (7) (2006): 1067- 
1080. 

Timoorei Sangani Ahmad. Calculation of the effective 
permeability and simulation of fluid flow in fracturated 
reservoirs, A dissertation submitted to the university of 
New South Wales (2005). 

Trochu, F. and Gauvin, R. Limitations of a Boundary-Fitted 
Finite Difference Method for the Simulation of the Resin, 
Transfer Molding Process, Journal of Reinforced Plastics 
and Composites., 11 (7) (1992): 772-786. 

Webb, Stephen W. Gas Transport in Porous Media, Part 1: 
Processes and Models, Edited by Clifford K. Ho, and 
Stephen W. Webb, Springer (2006): 55- 63. 


65 



